package biolabs.model.diffusion;

import java.util.List;

import biolabs.space.ContinuousPosition;
import biolabs.space.DiscretePosition;
import biolabs.space.GenericGrid;
import biolabs.space.Grid;

/**
 * Diffuser is an object that performs diffusion.
 * The diffusion coefficient and the diffusion mask are supplied
 * once at instance creation. The instance can then be used to
 * perform diffusion repeatedly.
 *
 * @param <P> position type of the gird
 * @param <C> continuous analog to {@code P}
 *
 * @author Tomas Mikula
 */
public class Diffuser<P extends DiscretePosition<P,C>, C extends ContinuousPosition<C,P>> {

	private final double coeff;
	private final Grid<Boolean, P, C> mask;
	private final List<P> stencil;

	/**
	 *
	 * @param mask positions which have {@code false} in {@code mask}
	 *   are omitted in calculation of neighbor concentrations.
	 * @param coeff diffusion coefficient. Valid range of values is
	 *   {@code (0; 1/(n+1)]} where {@code n} is the size of the
	 *   neighborhood of a grid element.
	 * @throws IllegalArgumentException if {@code coeff} has invalid value.
	 */
	public Diffuser(Grid<Boolean, P, C> mask, double coeff){
		this.mask = mask;
		this.coeff = coeff;
		this.stencil = mask.getSize().zero().getNeighbors();

		int stencilSize = stencil.size();
		if(coeff > 1.0/(stencilSize+1) || coeff <= 0)
			throw new IllegalArgumentException( coeff +
					" is invalid value for diffusion coefficient." +
					" Expected value from (0,"+1.0/(stencilSize+1)+"]");
	}

	/**
	 * @param conc grid of concentrations
	 * @throws IllegalArgumentException if sizes
	 *   of {@code conc} and {@code mask} differ.
	 */
	public void diffuse(Grid<Double, P, C> conc){

		P size = conc.getSize();
		if(!size.isEqualTo(mask.getSize()))
			throw new IllegalArgumentException("sizes of the concentration " +
					"grid and mask grid differ: "+size+" vs. "+mask.getSize());

		Grid<Double, P, C> buf = new GenericGrid<Double, P, C>(size.withChangedCoordinate(0, 2));
		calculateSubspace(1, buf, size.withChangedCoordinate(0, 0), conc, size.withChangedCoordinate(0, 0));
		int n = size.getCoordinate(0);
		for(int i=1; i<n; ++i){
			calculateSubspace(1, buf, size.withChangedCoordinate(0, i%2), conc, size.withChangedCoordinate(0, i));
			copySubspace(1, conc, size.withChangedCoordinate(0, i-1), buf, size.withChangedCoordinate(0, (i-1)%2));
		}
		copySubspace(1, conc, size.withChangedCoordinate(0, n-1), buf, size.withChangedCoordinate(0, (n-1)%2));
	}

	private void calculateSubspace(int nFixedCoords, Grid<Double, P, C> tgt, P tgtPos, Grid<Double, P, C> src, P srcPos){
		if(nFixedCoords == tgtPos.getDimension()){
			double c = src.at(srcPos);
			double diff = 0;
			for(P sp: stencil){
				P np = srcPos.plus(sp);
				if(src.inRange(np) && mask.at(np))
					diff += c - src.at(np);
			}
			tgt.setField(tgtPos, c - coeff*diff);
		} else {
			int n = tgt.getSize().getCoordinate(nFixedCoords);
			for(int i=0; i<n; ++i){
				srcPos = srcPos.withChangedCoordinate(nFixedCoords, i);
				tgtPos = tgtPos.withChangedCoordinate(nFixedCoords, i);
				calculateSubspace(nFixedCoords+1, tgt, tgtPos, src, srcPos);
			}
		}
	}

	private static
	<P extends DiscretePosition<P,C>, C extends ContinuousPosition<C,P>>
	void copySubspace(int nFixedCoords, Grid<Double, P, C> tgt, P tgtPos, Grid<Double, P, C> src, P srcPos){
		if(nFixedCoords == tgtPos.getDimension())
			tgt.setField(tgtPos, src.at(srcPos));
		else {
			int n = tgt.getSize().getCoordinate(nFixedCoords);
			for(int i=0; i<n; ++i){
				srcPos = srcPos.withChangedCoordinate(nFixedCoords, i);
				tgtPos = tgtPos.withChangedCoordinate(nFixedCoords, i);
				copySubspace(nFixedCoords+1, tgt, tgtPos, src, srcPos);
			}
		}
	}
}
